function [ blendSites, blendTraits, blendWindows ] = blendTraits( traits, sites, maxGap )
%Use traitsCondensed as input
%maxGap = 3;

maxTrait = max(max(traits));
numSites = size(sites, 1);
blendIdx = 0;

for curTrait = 1:maxTrait
    numTraitSites = 0;
    traitSites = [];
    for j = 1:numSites
       if (curTrait == traits(j, 1))
           numTraitSites = numTraitSites + 1;
           traitSites(numTraitSites, 1) = sites(j, 1);
       end
    end
    
    curTraitSite = -maxGap - 1;
    numCurSites = 0;
    for j = 1:numTraitSites
        % If we've hit a gap, we need to commit the current set to blend
        if (numCurSites ~= 0 && traitSites(j, 1) - curTraitSite > maxGap)
            avgSite = 0;
            maxSite = -1;
            minSite = -1;
            for k = 1:numCurSites
                avgSite = avgSite + curSites(k, 1);
                if (curSites(k, 1) < minSite || minSite == -1)
                    minSite = curSites(k, 1);
                end
                if (curSites(k, 1) > maxSite)
                    maxSite = curSites(k, 1);
                end
            end
            avgSite = round(avgSite / numCurSites);
            window = max(abs(maxSite - avgSite), abs(minSite - avgSite));
            blendIdx = blendIdx + 1;
            
            blendSites(blendIdx, 1) = avgSite;
            blendTraits(blendIdx, 1) = curTrait;
            blendWindows(blendIdx, 1) = window;
            numCurSites = 0;
        end
        
        numCurSites = numCurSites + 1;
        curTraitSite = traitSites(j, 1);
        curSites(numCurSites, 1) = traitSites(j, 1);
    end
    
    if (numCurSites ~= 0)
        avgSite = 0;
        for k = 1:numCurSites
          avgSite = avgSite + curSites(k, 1);
        end
        avgSite = round(avgSite / numCurSites);
        window = max(abs(curSites(numCurSites, 1) - avgSite), abs(curSites(1, 1) - avgSite));
        blendIdx = blendIdx + 1;
        blendSites(blendIdx, 1) = avgSite;
        blendTraits(blendIdx, 1) = curTrait;
        blendWindows(blendIdx, 1) = window;
        numCurSites = 0;
    end
    
end

end

